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1 Introduction 



Lord Kelvin was the first to suggest that knotlike configurations might be of fundamen- 
tal importance. In |]T] he proposed that atoms, which at the time were considered as 
elementary particles, could be viewed as knotted vortex tubes in ether. Subsequently he 
also conjectured [§], that (thin) vortex filaments that have the shape of torus knots f|, 
H should be stable. 

Kelvin's theory of vortex atoms has long ago subsided. However, at the time it 
was taken seriously, and it led to an extensive study and classification of knots. In 
particular the results obtained by Tait || remain a classic contribution to mathematical 
knot theory. Recently Kelvin's idea that different elementary particles can be identified 
with topologically distinct knots has been advanced in particular by Jehle ||. 

Since its inception, Kelvin's conjecture on the stability of torus knots has also at- 
tracted much interest. However, it seems that until very recently J7|, the issue has 
not been even properly addressed. This is due to the exceedingly complex nature of 
dynamical models that describe stable knots. 

In order to study the dynamics of knotlike configurations in a predictive manner, 
we need a Lagrangian field theory where knots can appear as solitons. But in such a 
theory the equations of motion are generically highly nonlinear, to the extent that an 
analytic approach is hopeless. Indeed, even the identification of a Lagrangian where 
stable knotlike solitons can be present, has until recently defied all attempts. With 
the recent, rapid advances in computers it appears that by combining topological and 
geometrical methods in quantum field theories with high performance computing, a first 
principles analysis of Kelvin's conjecture is finally a reality (7|. 

In addition of elementary particles (atoms), Kelvin was also interested in knotlike 
configurations in a hydrodynamical context. Today there are numerous additional sce- 
narios, where dynamical knots can be important: Presently it is commonly accepted that 
fundamental interactions are described by string theories 0, with different elementary 
particles corresponding to the vibrational excitations of a primary string. Even though 
connections between modern string theory and knot theory do exist, for example the 
Chern-Simons action is related both to conformal field theories and knot invariants 0, 
the possibility of an intimate relationship e.g. at some nonperturbative level remains to 
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be investigated. But in a number of other physical scenarios the potential relevance of 
knotlike structures has already been established ||, ||. For example in early universe 
cosmology cosmic strings are supposed to be responsible for structure formation. These 
strings are expected to decay e.g. via gravitational fluctuations, but the possibility that 
stable knotted strings might survive should have important consequences. Furthermore, 
in QCD one may expect that gluonic flux tubes that confine quarks in hadrons could 
become tangled, suggesting that in pure Yang-Mills theory closed knotted flux tubes 
appear as physical states. Stringy knotted vortices might also appear in a variety of 
condensed matter physics scenarios. For example in type-II superconductors magnetic 
fields are confined within the cores of vortex-like structures. Recent experiments with 
3 He-A superfluids have also revealed interesting vortex structures that can be described 
by theoretical methods which are adopted from cosmic string models. The study of knot- 
like configurations is also highly important for chemical compounds such as polymers. 
Finally, the investigation of knots is rapidly becoming an important part of molecu- 
lar biology, where entanglement of a DNA chain interferes with vital life processes of 
replication, transcription and recombination. 

Until now knotlike configurations have been mainly studied using non-dynamical, 
phenomenological models |[L0|| . One first introduces a one-dimensional, structureless 
string and distributes some charge along it. The configuration is then allowed to relax 
into an equilibrium shape. However, such an approach is not dynamical, there are in- 
evitably several knot-specific ad hoc parameters that need to be determined by various 
means. For example, the length of each topologically distinct knot is an independent de- 
gree of freedom. The details of the charge distribution also introduces some arbitrariness 
that can be parametrized. As a consequence of such parameters that are specific for in- 
dividual knot configurations, the ensuing models lack in their predictive power. There is 
a definite need for a first principles approach where different knots emerge dynamically, 
as solitons. One can then predict their properties in terms of fundamental quantities 
which are in no manner specific to a particular knotted structure but characteristics of 
the underlying physical environment. 

The literature on solitons is enormous, and there are several extensive reviews |JTTJ . 
Until now the activity has mainly concentrated on 1+1 dimensions with the notable 
exceptions of the 2+1 dimensional vortex and nonlinear cr-model solitons, and skyrmeons 
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and magnetic monopoles in 3+1 dimensions. These are all pointlike configurations, that 
can not be directly associated with knotted structures. 

When a pointlike two dimensional soliton is embedded in three dimensions, it be- 
comes a line vortex. Since the energy of a vortex is proportional to its length, for a finite 
energy the length must be finite. This is possible if the core forms a closed, knotted 
structure. In 1975 one of us [|1^] proposed that a torus-shaped closed vortex could be 
constructed in a definite dynamical model. The configuration suggested in [O 



is con- 



structed from a finite length line vortex, twisted once around its core before joining the 
ends. The twist ensures that the configuration is stable against shrinkage. As a knot it 
corresponds to the unknot which is the simplest possible knotlike configuration H], ||. 
By combining geometrical and topological methods with high performance computing, 
we recently verified [7[ the existence of such a soliton in the model proposed in [12 



Furthermore, we also found definite evidence for the existence of a soliton in the shape 
of a trefoil, which is the simplest possible torus knot. This is a strong indication, that 
the model proposed in |L2] realizes Kelvin's conjecture on the existence of stable torus 
knots. 

In this article we shall present a detailed investigation of the torus-shaped unknot 
soliton in the model proposed in fn| . In the next section we describe this model in detail, 
and in section 3. we specify it for torus-shaped solitons. In section 4. we reformulate 
the equations of motion so that a numerical construction of the unknot soliton can be 
performed, and in section 5. we present results of a model numerical simulation. 



2 A Hamiltonian with Localized Solitons 

In [0] we have proposed that certain relativistic field theories may admit solitons with 
the topology of a generic torus knot. In particular, we argued that such solitons can be 
described in terms of a three component vector field n(x) with unit length n • n = 1. 



The action is 12 



S = J d 4 x ^_L(,9 M n) 2 + ^(n-^nx3„n) 



(1) 
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where we choose n to be dimensionless so that e determines a length scale and g is a 
dimensionless coupling constant. This is the most general relativistically invariant action 
for n that involves second and fourth order derivative terms so that time derivates appear 



only quadratically. It can be viewed as a U(l) gauged version of the Skyrme model [] I 
However, the dynamical content of these two models are quite different. 

Since time derivatives in action (|IJ) appear only quadratically, it admits a canonical 
interpretation with static Hamiltonian 



H = E 2 + E 4 = Jd 3 x(^(d t n) 2 + ^(n-^nx^n) 



(2) 



The first term E 2 coincides with the nonlinear 0(3) cr-model action, known to admit 
stable solitons in two dimensions. But if we include the fourth order derivative term E4, 
stable finite energy solitons are also possible in three dimensions. This is suggested by 
the Derrick scaling argument: If we set 

x — > px (3) 

we find for our Hamiltonian 

H = E 2 + E 4 A P E 2 + -E A (4) 

so that stable finite energy solitons may exist when £4 is present. In particular, we 
conclude that such solutions obey the virial theorem 

E 2 = E 4 (5) 

A soliton described by n(x) is a localized configuration in R 3 , at spatial infinity 
|x| — > 00 the vector field n(x) approaches a constant vector no- Hence we compactify 
R 3 into a three dimensional sphere S 3 , and n(x) can be viewed as a mapping from 
the compactified R 3 ~ S 3 —>■ S 2 . Such mappings fall into nontrivial homotopy classes 
7r 3 (S 2 ) ~ Z that can be characterized by the Hopf invariant H], H. For this we introduce 
the closed two-form 

F = (dn A dn, n) (6) 
on the target S 2 . Since H 2 (S 3 ) = its preimage F* on the base S 3 is exact, 

F* = dA* (7) 
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and the Hopf invariant Qh coincides with the three dimensional Chern-Simons term, 

1 

4^ 



Qh = -^IFAA (8) 



R3 



The existence of stable solitons in (|2|) with a nontrivial Hopf invariant is then strongly 
suggested by the lower bound estimate 

H>c-\Q H \i (9) 



where c is a nonvanishing numerical constant |0| , 1TJ 



In the following we shall find it useful to introduce a unit four vector <& M (x) : R 3 ~ 
S 3 — > S 3 and define two complex variables 

Zx(x) = $ 1 (x) + ^ 2 (x) 
Z 2 (x) = $ 3 (x)+^ 4 (x) 

so that 

|Zi| 2 + |Z 2 | 2 = 1 
and describe the unit three vector n(x) by 

n k = Z ] a k Z (10) 

with a k the Pauli matrices. The two-form in becomes 

Fij = i(diZ^djZ - djZ^d.Z) (11) 

and modulo U(l) gauge transformations we have 

A^-i&diZ-di&Z) (12) 

Under a local U(l) rotation Z —> e n Z both F^ and n remain invariant while Ai suffers a 
gauge transformation Ai — > Ai+dtf. Hence Qh is a functional of n only, and substituting 
in (H) we get 

Qh = J d3x A d ®P A ( 13 ) 

which is the standard integral representation of the winding number for maps S 3 — > S 3 , 
invariant under arbitrary local variations $ M — > <J> M + 5 < I ) M . In particular, (§) is an integer 
valued invariant of the vector field n(x). 
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3 Toroidal Solitons 



In the following we shall be interested in the simplest possible knotlike soliton of (y), 
the torus-shaped unknot configuration ||12||, WM. In this case we can introduce further 



simplifications using the fact that the Hamiltonian (|2]) is invariant under both spatial 
SO (3) rotations and internal SO (3) rotations of the vector n(x). Since we expect a torus 
shaped soliton to exhibit rotation invariance around its toroidal symmetry axis, this 
global SO (3) x SO (3) invariance will be broken into a diagonal SO (2) G SO(2)xSO(2) 
invariance of combined spatial and internal rotations around the symmetry axis. Instead 
of minimizing the energy of the original Hamiltonian with respect to arbitrary variations 
of n(x), it is then sufficient to minimize the energy for the most general toroidal SO(2) 
G SO(2)xSO(2) invariant configuration n(x). 

We select the toroidal symmetry axis so that it coincides with the z- axis in R 3 , and 
align the asymptotic vector field n so that it points along the positive z-axis. In terms 
of cylindrical coordinates (r, z, ip) = (y, ifi), the most general Ansatz which is consistent 
with the toroidal SO (2) symmetry then separates the angle i)j that describes rotations 
around the z-axis, 

/ sinO(y) + kip) ■ sin6 l (y) \ 
n(y,V) = cos(^(y) + kif)) -sin0(y) (14) 
V cos^(y) ) 

Here k is an integer that counts the number of times the torus circles around the z-axis. 
Since (|14]) approaches the asymptotic constant vector no when |x| — > oo and no points 
along the positive z-axis, we conclude that 

0( y ) |y J^° mod 2vr (15) 

We substitute the Ansatz (|HD in (fj) and scale H appropriately by the dimensionless 
coupling constant g. In this way we find for the Hamiltonian 

H = G 2 -E 2 + E A = J drdz r ■ S(r, z) 



J drdzrl^(9 2 r +el) +sin 2 ^[^{^+^+^} + ^{(^A-^,) 2 + ^(^+^)} 



(16) 
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Here G is a dimensionfull coupling constant inversely proportional to length. Since there 
are no additional dimensionfull quantities, G specifies a length scale which determines 
both the size and the shape of the solitons. We note that H is invariant under z — > — z 
reflections. 

We introduce the following parametrization of the four- vector <&(x), 



/ cos 012 sin $ \ 
sin (pi2 sin $ 
COS 034 cos i? 
V sin cos $ J 



(17) 



From ([T0|), (fbl| ) we then find 



(p(y) + kip 



(P34 

2i) 



'12 



For the gauge field (|12[) 



and (|), (0) yields 



so that the Hopf invariant 



A = cos 2 1? (i0 3 4 + sin 2 $ i 

= <M = sin 2t? A (dcj)^ 
(M) becomes 



'12 



7Y^ 



Q H = — / s in2^ d-d A 



A 



>12 



(18) 

(19) 
(20) 

(21) 



Here the expected local S 1 x S 2 structure of our toroidal configuration is manifest. 

From (|16D we conclude that for a finite energy the angle 9 = 2d must vanish on the 
z-axis, 

9{r = 0,z) = mod 2n (22) 

Together with (|15D this implies that we can select 9 = on the entire boundary of the 
half-plane (r,z). In order to have a nontrivial Hopf invariant we need that at some 
interior points (r c , z c ) of the half plane we have 9(r c , z c ) = n mod 2n. For simplicity 
we assume that there is only one such point, and due to the z — > —z symmetry of flT5| ) 
we can select z c = 0. This point determines the center of the toroidal configuration, a 
circle of radius r c on the z = plane of R 3 which is parametrized by the angle ip. 
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In (1171) the preimage of 9 = corresponds to a circle which is parametrized by the 
angle ^34. Similarly, the preimage of 9 = n corresponds to a circle parametrized by the 
angle 0X2- Consequently we can refine (p"8|) by the further identifications 



<p(y) = 034 

kip = -0i2 (23) 

and for the Hopf invariant we obtain 

Q H = J sin 9 d9 A dtp A dip = - J sin 9 dB A dip (24) 

For a nontrivial Hopf invariant the angle *p{y) must then increase (decrease) by 2n (more 
generally by 2nn with n an integer) when we go once around the center r c on the (r, z) 
half plane. This means that the two circles parametrized by 0i2 and ^34 corresponding 
to the preimages of 9 = it and 9 = respectively, are linked (k ■ n times). 

Finally, if we assume that our toroidal soliton solution is an analytic function of the 
variables r and z, we conclude from ( |1~6| ) that its energy density vanishes at the center 
of the torus £(r c ,0) = 0, where 9 has a maximum value 9 = n. Since the minimum 
value 9 = occurs on the boundary of the (r, z) half plane, for an analytic soliton the 
energy density must also vanish there, and in particular along the z-axis. Hence S(r, z) 
must have a maximum value on a ring between the center r c and the boundary of the 
half plane. This means that the energy density of a toroidal soliton is concentrated in a 
tube-like neighborhood around its center. 

In the following we find it convenient to introduce stereographic coordinates U(y), 
V(y) on the (r, z) Riemann sphere, defined with respect to the north pole at 9 = 

V 

tp = — arctan(— ) (25) 



9 = 2 arctan VU 2 + V 2 (26) 
In these variables the Hamiltonian (|16|) becomes 

H = G 2 -E 2 + E A = J drdzr -£(r,z) 

r2 k 2 U 2 + V 2 



Jdrdzr\G\ — 2+ V 2)2 2 + 



(1 + U 2 + V 2 ) 2 r 2 (l + U 2 + v 2 y 
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(U r V z - U z V r ) 2 k 2 (UU r + VV r ) 2 + (UU Z + VV Z ) 2 \ 

(i + u 2 + v 2 y r 2 (i + u 2 + v 2 y j 1 ' 



and for the Hopf invariant we get 



_ k f dU A dV 
Qh ~ ~nJ (1 + U 2 + V 2 ) 2 (2 * 



4 Flow Equation 

The Euler-Lagrange equations for ( p^6| ) are highly nonlinear and an analytic solution 
appears to be impossible. It seems to us, that the only tools available are numerical. 
However, even a numerical integration is quite nontrivial, there are several complications 
that need to be resolved. One concerns the (expected) chaotic nature of the equations 
of motion: We do not expect a direct Newton's iteration to converge towards a stable 
configuration unless we succeed in constructing an initial configuration which is very 
close to the actual solution. Rather, we expect Newton's iteration to exhibit chaotic 
behavior. 

Due to lack of a proper initial configuration for Newton's iteration, we resort to 
alternative methods. For this we first formulate the problem at an abstract level, by 
considering a generic static energy functional E(q) with some variables q a . We introduce 
an auxiliary variable r, and extend the (stationary) Euler-Lagrange equations of E(q) 
to the following parabolic gradient flow equation 

^ = 

dr dq a 



Since 

dE /SE\' 
dr \6q a J 



(30) 



the energy decreases along the trajectories of (|2^). Furthermore, by squaring (|29|) and 
integrating from some initial value r = T to r — > oo we get 

T T 
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Hence r-bounded trajectories of ( p9|) flow towards a stable critical point of E(q). In 
particular, if we start at initial time r = T from an initial configuration q = q and 
follow a bounded trajectory of fl29|), in the r>T limit we eventually flow to a stable 
critical point of E(q) 

In (p?D the fields U(y), V(y) correspond to the variables q a , and by denoting Wi = U, 
W 2 = V the flow equation becomes 

dW a (y) _ 2 £E 2 ^4 , . 

5r " cW a (y) cW a (y) 1 ) 

If we introduce the scaling (H) we find for the scaled Hamiltonian the flow equation 

dW a (y) , ^2 ^4 , . 

«(Jr) " (pG) '^(y) 5W.(y) ( ' 3) 

Thus a flow towards a toroidal soliton with coupling constant G coincides with the flow 
towards a toroidal soliton with coupling constant p 2 G, provided we rescale the flow 
variable r into -r. This means that the soliton is essentially unique; It is sufficient 
to consider the flow towards a soliton with a definite value for G, since solitons with 
other values of G are obtained from this configuration by a simple scale transformation. 
Notice however, that in practice we perform a numerical integration of ( j3^) on a finite 
lattice, and there are finite size corrections to the simple scaling (j33|). 

An additional problem in a numerical approach is the selection of an optimal size 
for the lattice. If the lattice is too small in comparison to the scale that describes the 
soliton, the soliton may not fit into it. On the other hand, if the lattice is too large in 
comparison to the soliton, we may either miss the soliton entirely or use an unnecessarily 
large amount of computer time in constructing it. 

The coupling constant G in ( pTD is the sole dimensionfull quantity that appears in 
our equations. Consequently it determines how both the location r c of the center and 
the thickness of the soliton scale. Besides G, there are also dimensionless numerical 
factors that affect the overall size of the soliton. Since these numerical factors can only 
be obtained by actually solving the equations of motion, there are no a priori methods 
in selecting an optimal size of a lattice. 

In order to approach these problems, we have developed a simple renormalization 
procedure that allows us to select the initial configuration so that the location r c of its 
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center coincides with the location of the center for the actual soliton. For this we first 
observe that the soliton obeys the virial theorem (f|), in our present variables 



G 2 E 2 



E 4 




By demanding that this virial theorem is also obeyed during the flow ([32]), we promote 

the coupling constant G into a r-dependent variable G — > G(t), 



When we approach a soliton as r — > oo, the variable G(j) must then flow towards an 
asymptotic value G* which is the value of the coupling for the soliton, 



This renormalization procedure fixes the location r c of the center for the soliton, and 
allows us to choose the size of our lattice appropriately. 

The thickness of the final soliton is also determined by the asymptotic value G* of 
the coupling constant. However, since additional dimensionless numerical factors are 
also involved, the renormalization ([34] ) does not help us in selecting the thickness of 
the initial condition so that numerical convergence is secured. This poses a problem 
for which we at the moment lack a firm solution, besides experimenting with various 
different types of initial configurations. 

Finally, on a finite size lattice we also need to determine boundary conditions at the 
boundary of the lattice. On the (r, z) half plane the proper boundary condition for U (r, z) 
and V(r, z) is that both vanish on the boundary of the half plane. This corresponds 
to the compactification of the half plane into a Riemann sphere. However, we expect 
that on a finite lattice a simulation with such trivial boundary conditions leads to a flow 
towards the trivial configuration U(r, z) = V(r, z) = 0. In order to impose the boundary 
conditions properly in our numerical simulation, we have adapted an iterative process 
where we first specify the boundary conditions using the initial configuration £/o(y), 
Vo(y). At later values of r we then update these boundary conditions successively, by 
interpolating the iterated configurations from the interior of the lattice to its boundary. 
In this manner we expect that we eventually obtain boundary conditions that correspond 




G(t) 



■oo 



G* 



(34) 
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to those of an actual toroidal soliton. An alternative method that we have used to 
determine the boundary conditions, is to start the iteration of ([32|) from a sufficiently 
large initial lattice with boundary conditions determined by the initial configuration. By 
successively shrinking the size of the lattice and determining boundary conditions in the 
shrinked lattice using the appropriate restriction of the pertinent iterated configuration 
from the larger lattice, we expect to converge towards boundary conditions that coincide 
with those of an actual soliton. 



5 Numerical Simulation 

We shall now present an example of a numerical construction of an unknot soliton. For 
simplicity we restrict to a configuration with Hopf invariant Qh = 1 (i-e. k = n = 1), 
the other cases being treated similarly. 

For a numerical integration of (^), we need an initial configuration with Qh = 1 
that properly approximates the final soliton. For this we introduce the standard complex 
bilinear transformation 

« = 7^ < 35 > 

and view it as a map between the unit disk |£| < 1 and the Re(a) > half plane. By 
employing polar coordinates on the disk, we observe that the loci on the ty-plane are 
circles, that at least in a qualitative sense coincide with our a priori expectations for the 
angular fl25|) and radial (|26"D behavior of a soliton (in coordinates where the center of a 
soliton coincides with 6 = 0). 

Our soliton maps our compactified (r, z) half plane to a Riemann sphere. The vari- 
ables U (r, z) and V(r, z) are stereographic coordinates with respect to the north pole 

6 = 0, which topologically coincides with the entire boundary of the (r, z) half plane. 
As a consequence U and V determine a mapping from the (r, z) half plane to a disk, 
with the center of the disk at 8 = (the entire boundary of the (r, z) half plane), while 
the rim of the disk is an (infinitesimally) small circle around the center of the soliton at 
r = r c where we have 8 = n. This suggests that we start our construction of an initial 
configuration by first setting 

T -\- %Z -\- O 

W(r,z) = U(r,z) + iV(r,z) = — — (36) 

r + iz — r c 
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so that 

r — r, 



U(r,z) = 1 + (r c + p) ■ f , 



(r — r c ) 2 + z 2 

z 

r — r c )' z + z' z 



V(r,z) = -(r c + p). 7 —— —— (37) 



and the loci are circles on the U, V plane with 

r c + P, 2 (r c + p) 2 



(u-iy + (v + 



2z ' " 4z 2 



^ 1 1 r + P y + y2 l{r c + pY 



2 2 r — r c 4 (r — r c ) 2 

Here p > is a parameter, according to (j35l) the canonical value is p = r c but for the 
moment we leave it unspecified. 

For fl36f) we have 9(r c , 0) = ir in (^6|). But in addition we need 5 = on the entire 



boundary of the (r, 2) half plane. For this we improve ([$6|) into 

U(r,z) = Fu(r,z).(l + (r c + p)- {r ^ ^ ^ 

V(r,z) = -F v (r,z)-((r c + p). (r _ r * a + jg3 ) (38) 

where Fy, are appropriate positive valued functions, that vanish on the boundary of 
the (r, 2) half-plane; For a proper initial configuration both Fy and i*V vanish sufficiently 
rapidly when we move away from the center at r = r c , so that the initial configuration 
(in addition of the final soliton) can be fitted inside the lattice that we use in our 
simulation. We also need to specify these functions so that the Hopf invariant of the 
initial configuration coincides with the Hopf invariant of the desired soliton, in the 
present case Qh = 1- 

We have performed extensive numerical simulations for a soliton centered at r c = 6, 
using several different initial profiles Fjj and Fy, lattice structures and various alterna- 
tive methods for ensuring independence from boundary behavior. Our simulations are 
consistently converging towards a definite toruslike configuration. As an example we 
now describe results from a particular simulation, based on the initial configuration (|38| ) 
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with 



exp(— di(r — r c ) 2 — d,2Z 2 ^ 



Fu(r,z) = F v (r,z) = b — — (1 - tanWei(r - r c ) 2 + e 2 z 2 }) 

Ci + (r — r c ) z + ciz 1 v 

(39) 

where b, e 2 are parameters, fixed by minimizing the total energy ( p7[ ) with respect to 
these parameters. This minimization is subject to the supplementary condition, that the 
overall size and shape of the initial configuration should resemble as much as possible the 
fixed point configurations that we have found in our earlier numerical simulations. In 
particular, we use the angular (ip) profiles of our earlier simulations to select an optimal 
value for the parameter p in (|36|). 

In a simulation of the flow equation (|32l) we determine the length of each time step 
Ar_n = r n — r n _i adaptively, by demanding that the relative variation of total energy 
remains bounded, 

|*^| < 10- (40) 
' E(r n ) 1 " V ; 



We employ version 5.4 of the PDE2D finite element algorithm |15| using a sparse 
direct method with a 4 th order Hermitean polynomial finite element basis. We have 
chosen a finite element approach, since it computes a continuous piecewise polynomial 
approximation to the solution. In a problem of topological nature this should be a defi- 
nite advantage e.g. over a finite difference approach, where the solution is approximated 
only at discrete lattice nodes. In our example we have started from an initial square 
lattice with 0.001 < r < 14.0 and —6.0 < z < 6, divided into a finite element mesh with 
15.000 triangular elements. We have selected the triangulation so that it is more dense 
near r k and z ~ 0, where we expect to have strongest dependence on the initial 
boundary conditions. 

On the initial lattice we iterate the flow equation (|32"1) until we reach a fixed point 
configuration described in figures la, 2a, 3a. We then introduce a sublattice with 0.1 < 
r < 13.5 and —5.5 < z < 5.5, divided into a finite element mesh with 18.000 evenly 
distributed triangular elements. On this sublattice we determine both the initial and 
boundary conditions by restricting the previously obtained fixed point configuration to 
the sublattice. We then continue the iteration of the flow equation fl3"2] ) until we again 
reach a fixed point, described in Figures lb, 2b, 3b and 4-7. 

In this two-phase simulation each of the phases takes about 200 hours of CPU time, 
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on a Digital Alpha Server 8400 equipped with EV56/440 MHz processors and 4GB of 
CPU. In the second phase with a larger number of triangles, our computation consumes 
about 725MB of CPU. 

Our simulations reveal very strong convergence both for the coupling constant G(t) 
and for the total energy E(r) towards fixed point values; see figures 1-2. The length of 
the time step also increases exponentially when we approach the fixed point; see figure 3. 
Furthermore, we find very strong pointwise convergence for both of the angular functions 
9(r, z) and <p(r, z), and also for the density Qh(i~, z) of the Hopf invariant; see figures 4-6. 
The integrated Hopf invariant Qh is very stable under our entire simulation, essentially 
Qh ~ 1 and for the final configuration we have Qh = 0.99997..., when integrated over 
the final lattice. 

Finally, for the energy density £(r, z) in ( P7| ) we have also found very strong pointwise 
convergence except in the vicinity of the origin r ~ z ~ 0. In this region we experi- 
ence difficulties in fully eliminating the residual dependence of S(r,z) from our initial 
boundary condition along the z-axis. This slowness in pointwise convergence for very 
small r and z in £(r, z) can be related to our method: We use a finite element approach, 
where in each triangle of the mesh we approximate the solution by a 4 th order Hermitean 
polynomial in r and z. Since £(r, z) involves r~ 2 terms, it is then conceivable that if the 
derivatives of U(r, z) and V(r, z) in ( J?7D do not vanish sufficiently rapidly for small r, 
the pointwise convergence of £ (r, z) may become slow in this region. 

However, if we include the additional measure factor r in fl27D , it tames the r~ 2 
behavior to the extent that the integrand of the total energy does exhibit very strong 
pointwise convergence. This is also the reason, why the difficulties we experience in 
eliminating the boundary dependence from £(r, z) as r, z — > do not interfere with the 
overall convergence of the energy for our solution. 

In figure 7 we present the energy density S(r, z) of our final configuration for r > 0.3, 
which is the region where we have confidence in the pointwise convergence of £(r, z) in 
the present example. The result is consistent with our qualitative picture: The energy 
density has a definite tubelike shape, essentially vanishing inside a tubular region around 
the center at r = r c and then attaining a maximum value in a surrounding collar-like 
region. At z = the maximum value of the energy density is obtained for a relatively 
small value of r = r max , and by combining results from a number of different simulations 
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we have arrived at an estimate that r max < 0.9. For smaller values of r the energy density 
£(r, z = 0) then appears to decrease very rapidly, but due to the ensuing difficulties with 
small r convergence we can not verify that on the z-axis we actually have £{r — 0, z) ~ 0. 
For this one needs to perform a high precision three dimensional simulation, such a 
simulation is now in progress and we plan to report on it in a future publication. 

6 Conclusions 

In conclusion, we have analyzed in detail the torus shaped unknot soliton using a high 
precision simulation. We have found very strong convergence towards a critical point 
of the energy, suggesting that we have indeed found a very good aproximation to an 
actual soliton. Our results indicate, that a combination of geometrical and topological 
tools in quantum field theory with high performance computing is finally making a first 
principles investigation of knotlike solitons realistic. In particular, our results indicate 
that Kelvin's conjecture on the existence of stable torus knots can be realized in the 
model that we have studied. 
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Figure Caption 



Figure la: Flow of coupling constant G(t) towards its fixed point value G*, for a 
soliton with r c = 6. 

Figure lb: Flow of coupling constant G(r) towards its fixed point value G*, starting 
from the fixed point configuration described in Fig. la on a truncated sublattice, with 
boundary conditions also determined by the fixed point configuration in Fig. la. The 
slight difference in the fixed point values G* after the first and second phase is partly due 
to a smaller and relatively more accurate lattice employed in the second run, and partly 
due to an improvement in the boundary conditions. Due to improved convergence there 
is also a O(10 3 ) increase in the final time step scales. 

Figure 2a: Flow of total energy as a function of r towards a fixed point. Notice 
that the energy decreases as a function of r, as it should for a bounded trajectory. 

Figure 2b: Flow of total energy as a function of r towards a fixed point, starting 
from the fixed point configuration described in Fig. 2a on a truncated sublattice, with 
boundary conditions determined by the fixed point configuration in Fig. 2a. Notice that 
there is a O(10 3 ) improvement in the time step scales. 

Figure 3a: Length of time step as a function of the number of iterations. The length 
of a time step is determined by the bound fl40f) . The length should increase without a 
limit when a fixed point is approached. 

Figure 3b: Length of time step for the second phase, starting from the configuration 
described in Figure 3b. The relative increase by a factor of O(10 3 ) from Fig. 3a can be 
attributed to improvements in boundary behavior and lattice accuracy. 

Figure 4: The angle 9(r, z) for a fixed point configuration at the end of the second 
phase. Plotted on a sublattice with 0.3 < r < 12 and —4.5 < z < 4.5. At the center 
r c = 6 we have 9 = tt, and outside of the soliton 9 ~ 0. 

Figure 5: The angle <p(r, z) for a fixed point configuration at the end of the second 
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phase. Plotted on a sublattice with 0.3 < r < 10 and —4.5 < z < 4.5. When we go 
once around the center at r c = 6, ip jumps by 2ir. Notice the appearance of a tubular 
structure in the middle, corresponding to the interior of the soliton. 

Figure 6: The density Qh{^, z) of the Hopf invariant, for the fixed point configura- 
tion at the end of the second phase. The tubular structure of a soliton is clearly visible. 
The integrated value of the Hopf invariant is very stable under our entire iteration, and 
for the final configuration we have Qh = 0.999997... 

Figure 7: The energy density S(r, z) for the fixed point configuration at the end of 
the second phase, plotted on a sublattice with 0.3 < r < 10 and —4.5 < z < 4.5. The 
tubular structure is clearly visible, but near r z we still have residual dependence 
on the initial boundary on the z-axis. This is understandable, since even a small error 
in the boundary condition can eventually give rise to a large local tension. 
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Figure lb 
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